Neuron dynamics in the presence of 1// noise 



Cameron Sobiejj Arif Babul|j and Rogerio de Sousc[^ 
Department of Physics and Astronomy, University of Victoria, Victoria, B.C., V8W 3P6, Canada 

(Dated: April 14, 2011) 

Interest in understanding the interplay between noise and the response of a non-linear device 
cuts across disciplinary boundaries. It is as relevant for unmasking the dynamics of neurons in 
noisy environments as it is for designing reliable nanoscale logic circuit elements and sensors. Most 
studies of noise in non-linear devices are limited to either time-correlated noise with a Lorentzian 
spectrum (of which the white noise is a limiting case) or just white noise. We use analytical theory 
and numerical simulations to study the impact of the more ubiquitous "natural" noise with a 1// 
frequency spectrum. Specifically, we study the impact of the 1/f noise on a leaky integrate and fire 
model of a neuron. The impact of noise is considered on two quantities of interest to neuron function: 
The spike count Fano factor and the speed of neuron response to a small step-like stimulus. For the 
perfect (non-leaky) integrate and fire model, we show that the Fano factor can be expressed as an 
integral over noise spectrum weighted by a (low pass) filter function given by J-{t, f) = sinc^(7r/t). 
This result elucidates the connection between low frequency noise and disorder in neuron dynamics. 
Under 1/f noise, spike dynamics lacks a characteristic correlation time, inducing the leaky and 
non-leaky models to exhibit non-ergodic behavior and Fano factor increasing logarithmically as a 
function of time. We compare our results to experimental data of single neurons in vivo [M.C. Teich, 
C. Heneghan, S.B. Lowen, T. Ozaki, and E. Kaplan, Journal of the Optical Society of America A 
14, 529 (1997)], and show how the 1/f noise model provides much better agreement than the usual 
approximations based on Lorentzian noise. The low frequency noise, however, complicates the case 
for information coding scheme based on interspike intervals by introducing variability in the neuron 
response time. On a positive note, the neuron response time to a step stimulus is, remarkably, 
nearly optimal in the presence of 1// noise. An explanation of this effect elucidates how the brain 
can take advantage of noise to prime a subset of the neurons to respond almost instantly to sudden 
stimuli. 

PACS numbers: 87.19. L-; 87.19.1c; 05.40.-a 



I. INTRODUCTION 

One of the major puzzles of neuroscience is how neu- 
rons can store, process, and compute despite the fact 
that the brain is extremely noisy [I] . Understanding the 
evolved mechanisms and the associated non-linear dy- 
namics that allow the neurons to function in - and even 
exploit - a noisy environment is an essential step to- 
wards gaining insight into the information transmission 
and communication networks in the brain. Such stud- 
ies also have important implications beyond the domain 
of biophysics and neuroscience. Noise provides a critical 
barrier for the development of sensitive electronic and 
mechanical devices, particularly at the nanoscale [HE]. 
Increasingly, researchers are focusing on exploring inno- 
vative, non-traditional device design and control strate- 
gies that exploit the ambient noise [3HTT| . In this regard, 
there are clear advantages to understanding how nature 
has managed to harness noise in a setting whose primary 
(apparent) function is to manage information. 

It is generally accepted that neurons communicate 
with each other using sharp electric pulses referred to 
as action potentials or spikes. Each neuron is connected 
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to several other neurons, and will only generate a spike 
output when the integrated input from other neurons 
exceeds a certain threshold [T^l - fT^ . A startling discov- 
ery that (under certain circumstances) neurons can spike 
more regularly when stimulated by noise [15H17) led to 
assertions that noise is inherent to neuron function. Sev- 
eral subsequent experimental and theoretical studies were 
aimed at elucidating the functionality of neural noise |18l - 
[25] . In the first instance, noise - as expected - introduces 
a variability in the interspike intervals and degrades the 
information capacity of the spike trains, with low con- 
trast signals being most affected [55]. At the same time, 
these studies also found that stochastic resonance pro- 
vides a mechanism for neurons to take advantage of their 
own noise. In stochastic resonance, the addition of an 
appropriate amount of noise in a non-linear system can 
induce regularity by sensitizing subthreshold excitations, 
thus providing the extra energy for them to reach thresh- 
old [TSJ Uni [E] and enabling their detection. Addition- 
ally, Brunei et al. [21] and Svirskis [5S| have shown that 
a model neuron, when subjected to low frequency noise, 
is able to respond faster to a sudden excitation than in 
the absence of noise. For an animal living in a natural en- 
vironment, the ability to react quickly to sudden threats 
can mean the difference between life and death. 

All of the studies to date that have considered the im- 
pact of low frequency noise on neurons tend to model 
noise characterized by a single Lorentzian power spec- 
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trum. Natural noise, however, has an ubiquitous 1// 
frequency dependence [31 [26H30] . From current carrying 
electronic devices and geophysical time series to biolog- 
ical systems, the 1// power spectrum is everywhere. In 
biological settings, human hearing and speech [57], the 
response of biological photoreceptors to large intensity 
variation of visual image streams in nature |28j , the stride 
intervals time series of normal human gait [29| , intrinsic 
noise in neuronal membranes due to stochastic opening 
and closing of the various ion channels [30], etc. all ex- 
hibit 1/f behavior. Guerra et al. [7] demonstrated how 
stochastic resonance induced by 1/f noise can increase 
the sensitivity of nanomechanical resonators, allowing for 
the possibility of fashioning them into noisy but robust 
nanoscale computation devices. Neither the response nor 
the details of the underlying non-equilibrium behavior 
that makes neurons robust to natural {1/f) noise is well 
understood. 

In this paper, we present a comprehensive study of how 
neuron dynamics is affected by an arbitrary noise spec- 
tral density, and which sectors of the spectra are respon- 
sible for the beneficial functions that noise can provide. 
Specifically, we explore whether neuron response under 
1/f noise is significantly different from the one found in 
the presence of a simple Lorentzian spectra. 

The relevance of low frequency noise, implied by a 1// 
spectrum, to spike dynamics at the single neuron level 
is especially evident, as we shall argue, in the direct ex- 
perimental measurements of the spike count Fano factor 
by Teich et al. (HU 112] (Fano factor is the ratio between 
the variance and the mean spike number during a given 
observation time) . These authors demonstrated that the 
Fano factor of single neurons in the visual systems of 
cats and insects increases monotonically as a function of 
time. This is in dramatic contrast to the simple Poisson 
model (white noise), which leads to Fano factor equal to 
one at all times. The monotonic rise is also incompati- 
ble with models based on Lorentzian noise because the 
resulting Fano factor saturates at times longer than the 
inverse Lorentzian half width [3TJ [33] . We show that the 
characteristic non-ergodicity of 1/f noise explains why 
the Fano factor never saturates in single neuron experi- 
ments. Moreover, the rate at which the Fano factor grows 
as a function of time is different for 1/f and Lorentzian. 

In addition, we consider the effect of 1// noise on the 
reaction time of a neuron in response to a sudden stim- 
ulus. We demonstrate that 1/f noise is nearly optimal 
for speeding up neuron response. We provide an expla- 
nation for this effect that sheds light on the mechanism 
of neuron adaptation to their noisy environment. 

This paper is organized as follows. Section |lT| describes 
our model for the neuron: the leaky integrate and fire 
(LIF) model of the neuron, and explains how we intro- 
duce 1/f noise and other spectral densities in this model. 
Section [TTT| describes an analytical theory and a set of nu- 
merical simulations of the neuron Fano factor as a func- 
tion of time and compares our result to experime nt |31| . 
Particularly notable is our general expression Eq. (|13[) re- 



lating the Fano factor to an integral over low frequency 
noise weighted by an appropriate filter function. Sec- 
tion |IV| addresses the question of how noise can provide 
a mechanism for neurons to respond faster to a sudden 
stimulus. Section |V] provides our concluding remarks. 



II. THE LIF MODEL 

From a biophysical perspective, the classical Hodgkin- 
Huxley model [13,, £md its contemporary variants repre- 
sent the most realistic mathematical description of elec- 
trical response of a single neuron. Due to their intrinsic 
complexity, however, such models render the theoreti- 
cal and computational analysis of neuronal and neural 
network dynamics exceedingly difficult. For this reason, 
most studies to date tend to reference the simpler spik- 
ing neuron models, of which the leaky integrate and fire 
(LIF) model [12, JJ:^ that wc adopt is one. The LIF model 
represents each neuron by an electrical circuit; when 
appropriate circuit parameters and features are chosen, 
the LIF model can reproduce quite similar dynamics to 
the one described by the more complex Hodgkin-Huxley 
model [Tj. 

The LIF model consists of a capacitor C in parallel 
with a resistor i?; an injected continuous current I{t) 
models the spike input from a large number of neighbor- 
ing neurons. The neuron (or capacitor) voltage V{t) is 
given by the circuit equation. 



(1) 



A spike is generated whenever the voltage across the ca- 
pacitor reaches a certain threshold Vth! after the spike 
is emitted, the neuron is reset to a zero voltage state. 
Note how the threshold rule for spike generation intro- 
duces non-linearity in the LIF model: Consider two in- 
put currents Ii{t) and l2{t)] a neuron subject to input 
[Ii{t) -\- l2{t)] will generally reach threshold faster than 
a neuron that is subject to either Ii{t) or l2{t) only. 
Hence, the sum of outputs obtained from /i(t) and l2{t) 
applied separately is different from the output obtained 
from {Ii{t) + l2{t)]- Also, the resistance R plays an im- 
portant role in the model: it allows charge to leak out, 
thus negating inputs received in the distant past. Neu- 
rons have the property that inputs long past have less 
effect than recent inputs; a sufficient amount of input 
must happen sufficiently rapidly for the neuron to fire. 

The version of the LIF model that we are using has 
an additional feature that makes it more realistic: The 
introduction of a refractory time period r^, which models 
the physical reset time for a neuron after emitting a spike. 
This prevents the neuron from receiving input for a time 
Tr after spiking. Our choice for these circuit parameters 
are given in Table [l] 

The presence of a leak and a refractory period makes 
the LIF neuron extremely hard to treat analytically 03] ■ 
As a result, many theoretical studies have focused on the 
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R = Qo and Tr = limit, the so called perfect (non-leaky) 
integrate and fire model [2TJ [34]. This latter model is 
much easier to analyze but as we will demonstrate, the 
reduced complexity also leads to significantly different 
dynamics. 



A. Introducing noise in the LIF model 

We considered the LIF model subject to a noisy input 
current of the form 



i{t)^e{io + iiv{t)).. 



(2) 



where /q is a (constant) bias current, Ii is the noise am- 
plitude, and d{x) is the Heaviside step function: 9{x) = 1 
for a; > and 9{x) — for a; < 0; this ensures the current 
input represents the sum of spikes from a large number 
of connected neurons. The time series r]{t) is a Gaussian 
stochastic process, with variance equal to one and power 
spectra given by 
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with the brackets ( ... ) denoting ensemble averages 
over a large number of time series r]{t). Appendix A de- 
scribes the method used to generate individual time series 
for any given noise spectral density S{f). We considered 
a number of different noise densities, including the family 
of power law spectral densities 



1 



(4) 



where a is an exponent (a — 1 corresponds to 1// noise). 
The normalization constant Aq, is set by the condition for 
the variance to be one, 2tt J dfS{f) = {rf) = 1 [In the 
case of a = 1 Eq. Q is valid for 7,ni„ < / < 7max; 7min 
is a lower cut-off for which Sa{f) saturates, and 7max is 
an upper cut-off for which Sa{f) goes to zero faster than 
1/ P - see below]. 

Another important class of noise spectral density arises 
when the environment fluctuates with a single character- 
istic time Tc'- 



S{f) = 



2^2 /2 
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with 7 — \/[2nTc). This is a Lorentzian power spectrum 
and it implies S{t) = {r]{t)r]{0)) = e"!*!/^"^. Many au- 
thors refer to Tc as the "correlation time", and to the 
Lorentzian spectra as "time correlated noise". In the 
limit that 7 goes to infinite (tc — >■ 0), S{f) w {2tt'^j)~^ is 
approximately constant for all / <§; 7. Hence 7 — >■ 00 is 
the "white noise" limit. Another important limit occurs 
when 7 — 7^ (tc — >■ 00): In this case S{f) — > ^5{f)i 
signaling a "static" limit. 

It is useful to recall the basic physical picture for the 
origin of 1// noise. It emerges from the combination of a 



large number of Lorentzian fluctuators with an exponen- 
tially wide distribution of characteristic rates 7 [21 [3]. For 
example, assume 7 — 7maxe~^, with A a random variable 
that represents a distribution of activation energies. As- 
suming A is uniformly distributed in the interval [0, Amax] 
we get 



~S{f) 



dX I 7 

1 /-Tmax ^ 



27r2A, 



(6) 



arctan (j^^^ ~ ^-^ctan (^-^^ 



27r2A, 
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When 7„ii„ <C |/| <C 7maxj we may approximate 
arctan(7,„ax//) ~ 7r/2 and arctan(7,„in//) « 0; this 
leads to 
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where we used the fact that Amax = ln(7max/7min)- 
Hence, overall the resultant noise is well described by 




< I/I < 7min 
Tmin ^ I y I ^ Tma:5< 
Tmax 

< I/I < 00 



(8) 



with constant Ai = [47rln (7inax/7min)] ^ 

From Eq. ^ we see that the distribution of Lorentzian 
linewidths 7 is given by ^(7) — l/\dj/d\\ = 1/7. This 
is the reason why the spectrum acquires the 1// depen- 
dence. We may generalize this distribution to ^(7) — 
1/7", with a a dimensionless exponent; carrying though 
a similar derivation as in Eq. ^ leads to S{f) oc 1//". 
This shows that deviations of the I/7 distribution will 
reflect directly into a a 7^ 1 exponent for the noise spec- 
trum. 

Usually, 7,nin is exponentially small, and the experi- 
mental observation time window T is smaller than 7";^. 
In this case the low frequency cut-off will be instead set 
by 7inin — T^^. As the observation time T increases, 
more low frequency fluctuators will play a role; as a re- 
sult, I// noise has no characteristic time scale, and dis- 
plays non-ergodic behavior (time averages of observables 
are non-convergent, and can not be equivalent to ensem- 
ble averages). Below we discuss how the non-ergodic 
property leads to an increasing neuron Fano factor as 
a function of time. 



B. White versus 1// noise in the superthreshold 
regime: Bursting phenomena 

To compare LIF dynamics under the effect of white 
and 1// spectra, we considered noise in the superthresh- 
old regime (/q > Vth/i? = 4.28 x 10" A). We assumed 
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TABLE I. Circuit parameters used in our LIF model. The 
parameters are similar to the ones used to describe neurons 
in the cat's visual cortex (Chapter 14 of W). 



Resistance 


R = 38.3 Mfi 


Capacitance 


C = 0.207 nF 


Circuit time constant 


RC = 7.93 ms 


Threshold voltage 


Vth = 16.4 mV 


Refractory period 


Tr = 2.68 ms 



Iq = 4.3 X 10~^° A, slightly above threshold, ensuring 
that without noise the neuron will spike every 46 ms. 
For the cases with noise, we used Ii = O.l/p. The al- 
gorithm of Appendix A was used to generate a current 
input I{t), and Eq. ([I]) is integrated using the Runge- 
Kutta method. Figures |l(a) and |l(b) shows the neuron 



voltage as a function of time, for a particular time se- 
ries (the observation time window was T = 2 s). Under 
white noise, the spiking remains quite regular over time, 
because I(t) varies rapidly and most of its fluctuating 
components are filtered out. On the other hand, 1/f 
noise shows a combination of long periods of inactivity, 
with the voltage taking a long time to reach threshold, 
together with periods of spike bursting where the voltage 
reaches threshold on a much shorter time scale. This is a 
result of the fact that under 1/f noise, the current tends 
to get "stuck" at either small or large values. 

Figure [2] shows the interspike time interval histogram 
(ISI) for a 100,000 ensemble of time series with the same 
parameters considered in Figs.[T(a) and 1(b) Here we see 
that both types of noise can cause a notable decrease in 
the mean interspike time interval (compare to the noise- 
less case of a constant bias current). However, 1/f noise 
leads to a far more dramatic shift. A large portion of 
the interspike time histogram lies significantly below the 
noiseless interval, and a long tail is observed at large in- 
terspike times. This behavior is characteristic of bursting 
(See e.g. Chapter 16 of [1]); several spikes occur in rapid 
succession, followed by a longer period without spike ac- 
tivity. 

Figure [2] also shows the ISI for two Lorentzian noise 
spectra, with 7 = Hz and 7 = 1/T = 0.5 Hz, where 
T = 2 s was the simulation time window for each time 
series. Lorentzian noise with such low frequency leads 
to an ISI that is nearly as broad as the 1/f noise case; 
however, the Lorentzian noise cases do not display the 
long time tail characteristic of 1// noise. The simulations 
for 7 = Hz can be compared to an exact analytical 
result obtained using the methods of Refs. [2TJ |34] (See 
Appendix [b]) . Note how the 7 = Hz simulations are in 
excellent agreement with the exact result. 



III. NEURON FANO FACTOR UNDER 1/f 
NOISE 



Fano factor in the perfect integrate and fire 
model: Analytical results 



The Fano factor is defined by[T] 

{[ANit)] 



{N(t)) 



(9) 



where the random variable N{t) is the number of spikes 
generated from t' = to f = t, and AN{t) = N{t) - 
{N{t)). Hence, the Fano factor measures the amount 
of uncertainty in the spike train at a given time t. A 
noiseless spike train with identical interspike time inter- 
vals yields F{t) — 0. In contrast, consider the case that 
the spike events are uniformly distributed in the interval 
[0,t]: In this case the probability for a spike event to hap- 
pen during a time interval [f , t' + At] is independent of t' 
and given by iiAt, where fi is the mean firing rate. Then 
the resulting N{t) is a Poisson random variable uncorre- 
lated in time (white noise) with {N{t)) = {{AN)^) = ^i, 
leading to F{t) = 1 at all times. In the presence of low 
frequency noise, F{t) is known to become larger than one 

[211 [33]. 

The observation of a Fano factor much larger than one 
rules out the simple Poisson model, and suggests the pres- 
ence of long time correlations in the data [U [3T] . Mid- 
dleton et al. [21] derived an analytic expression for the 
Fano factor of the perfect (i? = 00) integrate and fire 
model subject to Lorentzian noise [Eq ([5])]. At times 
much longer than the average interspike interval their 
result becomes 
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where ((A/ )) = is the current variance, (/) — Iq 
is the bias current, Tc is the correlation time of the 
Lorentzian noise, C is the capacitor's voltage, and Vth 
is the threshold voltage. We emphasize that Eq. ( 10 ) as- 
sumes an input current /(t) = Iq + Ii-qit), i.e., it neglects 
the step function used in our numerical computations 
[compare to Eq. 

Here we generalize this result to an arbitrary noise 
spectral density. An RC circuit with no leakage {R = 00) 
will lead to a capacitor voltage V that always increases 
with increasing time. In the case of the neuron, V is re- 
set to zero when it reaches the threshold Vth- In other 
words, the voltage is decreased by Vth each time the neu- 
ron spikes. An equivalent way to treat this reset process 
is to instead increase the threshold by an additional Vth 
each time the neuron spikes; this allows us to count the 
number of spikes at a given time t by simply dividing 
the monotonically increasing V by Vth- Therefore the 
random variable N{t) is well approximated by 
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(11) 
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FIG. 1. (color online). Neuron voltage as a function of time (black curve), for (a) white noise and (b) 1// noise. The vertical 
lines (red) denote spiking events. Note how 1// noise leads to long time intervals with no spike generation, followed by intervals 
with spike bursting. 
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online) Interspike interval histogram 
of white and 1/f noise considered in 
respectively. The noise amplitude was 
h — O.l/o, with bias current lo slightly above threshold. Also 
shown is the case of no noise (/i = 0) , and of Lorentzian noise 
with 7 = Hz and 7 = l/T = 0.5Hz, where T = 2 s was the 
simulation time window for each time series. Note how 1/f 
has much broader ISI, with a long time tail extending well 
beyond its mean interspike time. 



Thiis approximation is valid at long times, t 3> CVth/Io, 
i.e. times much longer than the mean interspike interval 
so that N{t) can be represented by a real number instead 
of an integer. 

The simplicity of the perfect integrate and fire model 
lies in the fact that we do not need to consider the thresh- 
old barrier explicitly. This property relies heavily on the 
fact that the voltage of an RC circuit never decreases 
when R = 00. It is worthwhile to show how the approx- 



imation Eq. (11) fails in the presence of any amount of 
leakage. 

When R < 00, the mean voltage according to Eq. ([T]) 
is given by RIo{l — e~*^^'~'). Hence when t 00 the 
mean voltage saturates at RIq, and the approximation 
of Eq. ([n]) would give {N{t 00)) = RIo/Vth < 00. 
This is clearly an unphysical result: If the neuron spikes 
at least once, it will spike an infinite number of times 
when t — >■ 00, and {N{t — >■ 00)) must be either or 00. 
This unphysical saturation of {N{t)) shows why we must 
include the threshold barrier explicitly when dealing with 
the leaky model; it also shows why it is so difhcult to treat 
the leaky model analytically. 

In the absence of leakage, we have {N{t)) = 
t lo / [CVth) , i.e. at long times the number of spikes in- 
creases indefinitely with increasing time. The Fano factor 
can be calculated explicitly plugging Eqs. ^ and (11) 
into Eq. 
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(12) 



In the last step we changed the variables to T = {t'+t") /2 
and T = {t' — t"), and used the symmetry S{t) = S{—t). 
Eq. (12) provides an explicit relationship between the 



Fano factor and an arbitrary time correlation function. 
For example, upon inserting S{t) = c^l'^l/'^" we recover 
Eq. (fTol) exactly 



A simpler expression can be derived by inserting S{t) 
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2tt J dfe-'^^'f^Sif) into Eq. ^ 
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dfSif)T{tJ), 



where we defined the filter function by 

=sinc2(7r/t), 



(13) 



(14) 



where sinc(x) = sin (x) jx is the "unnorniaUzed sine func- 
tion" . Hence the Fano factor can be expressed as an inte- 
gral over low frequency noise; the filter function J^(t, /) 
dictates how much noise is "allowed in" at each given 
time t. Inspecting Eq. ( 14 ) shows that it acts as a low 
pass filter with bandwidth « Xjiixt). In the limit i — > oo, 
we have /t 3> 1 for all frequencies, and T « 5{f)/t. This 
leads to a useful result 



F{t ^ oo) = 
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(15) 



Therefore, the saturation of the Fano factor (or lack 
thereof) at long times is directly proportional to the 
amount of zero frequency noise. 



We can use Eq. ( 12 ) to find the neuron Fano factor 
subject to 1// noise. From Eq. ^ we know that the 
time correlation function for 1// noise can be written as 
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Hence the Fano factor for 1// noise can be written as a 
weighted average of Lorentzian Fano factors, 
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where Ce = 0.5772 is the Euler-Mascheroni constant. 
The latter approximation is valid for ^^^^ <^ t ^ — 
T, where T is the experimental time window. Hence 
Fi/f{t) increases monotonically with increasing time, un- 
til it reaches a saturation value around ~ T / In (T). This 
saturation value, however, is artifact of the finite length T 
of the experiment; it diverges as T — >■ cxi, i.e., when longer 
data sets are acquired. This effect is a manifestation of 
the non-ergodicity of 1// noise. In stark contrast, for an 
experimental window T > t^, the Fano factor for the sim- 
ple Lorentzian noise saturates at FLor. — 2Tclf / (CVth/o)j 
and remains at this value regardless of whether T is 
increased further. Assuming that this behavior carries 
over to the Leaky model, it offers one way to determine 
whether the neuron noise is better described by 1// or 
Lorentzian spectrum. The dependence of F{t) on t is 
another potential discriminator, as we shall see below. 

Figure [3] illustrates the relevance of the filter function 
Eq. (|14|) in quantifying the amount of noise absorbed by 
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FIG. 3. (color online) Fano factor [Eq. Q] for the case R = oo 
(non-leaky or perfect) integrate and fire model as a function 
of time for white noise, Lorentzian noise with 7 = 0.16 Hz, 
and 1 // noise with diflferent observation time windows T. For 
each noise spectra, the noise ampUtude Ii was chosen so that 
F(t = 0.1 s) — 10 for all noise spectra. This ensures that 
at t = 0.1 s the neurons absorb the same amount of noise 
power, despite the fact that the noise spectra are quite dif- 
ferent. Nevertheless, for times t after 0.1 s, the Fano factor 
diflers considerably for different spectra. This occurs because 
neurons "integrate noise" over a ban dwidth ~ 1/ {nt) as de- 
scribed by the filter function Eq. (14 1. 



neurons at a given time t. Here we plot the Fano factor 
F{t) for several different noise spectra, but choose each 
noise power (proportional to If) so that the Fano factor 
at a particular time t — 0.1 s is identical [-F(0.1 s) — 10] 
for all noise spectra. Hence at this particular i = 0.1 s 
the amount of disorder on neuron response is the same 
even though we are describing neurons subject to very 
different dynamical environments. Nevertheless, at times 
after t = 0.1 s the Fano factor differs considerably for 
different environments. This is a direct consequence of 
the fact that neurons integrate noise over a bandwidth 
~ l/^nt); hence, as t increases, a neuron absorbs noise 
over an increasingly narrow frequency range. Note also 
how F{t) for 1/f noise depends on the total observation 
time window T, and how F(t) is sensitive to different 
noise spectra before t — 0.1 s. 



B. Fano factor in the perfect integrate and fire 
model: Numerical results 



Figure |4] shows the Fano factor, as a function of time, 
for /o = /i = 2 X 10-1° A and T = 10^ s for white noise, 
and for Lorentzian noise with 7 = 1 Hz [t^ = l/(27r) s] 
and 7 = Hz (tc = 00). We also plot the analytical ex- 
pression Eq. ( 10 1 for the 7 = 1 Hz Lorentzian Fano fac- 



tor; as expected, the analytical expression shows slightly 
higher disorder [larger F(t)] than our Lorentzian noise 
simulation, because the latter only includes positive in- 
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FIG. 4. (color online) Fano factor [Eq. Q] for the R ^ oo 
(non-leaky or perfect) integrate and fire model as a function 
of time for white noise, Lorentzian noise with 7 = 1 Hz, and 
7 = Hz, the static case. We assumed Iq = h = 2 x 10~^" A 
and other parameters as in Table [l] Also shown is a compar- 
ison between the numerical simulation and the analytical ex- 
pression Eq. ( |10[ ) for Lorentzian noise. For white noise, F{t) 
tends to a quite small value at long times; for Lorentzian 
noise, F{t) plateaus at times t ~ (2^7)"^^. 



put currents (they are in close agreement when Ii ^ Iq). 
For white noise, the Fano factor tends to a small value 



at long times, in accordance with Eq. ( 15 ) that gives 
F{t 00) cx l/7max (7max = 10^ Hz is the upper fre- 
quency cut-ofF of our white noise spectrum). For the 
7 = Lorentzian, we have F(t) cx t in accordance with 
the hmit Tc -J> oo of Eq. ([lO| . 

Figure [5] shows the neuron Fano factor subject to 
1/f noise, using two different simulation time windows: 



T = 10^ s and T = lO'^ s. 



We also show the corre- 



sponding analytic results using Eq. ( 17) with 7,jii 
and 7max = N/T {N is the number of frequency inter- 
vals used in our simulation, see Appendix [A|) . Similar to 
the Lorentzian case, the analytic expressions for F{t) are 
larger than the numerical results because the former does 
not take into account the step function in Eq. ([2|. 

As expected, we find that the qualitative behavior of 
1/f noise is markedly different from Lorentzian noise. 
While for Lorentzian noise F{t) increases linearly with t 
until it reaches an asymptotic maximum at t « Tc, for 
1/f noise F{t) increases logarithmically [oc —tln(t/T)] 
and only reaches a slight saturation when t ^ T, the 
maximum possible value of time. 



C. Fano factor in the leaky integrate and fire model 

We now describe the impact of leakage on the Fano 
factor. We simulated the LIE model under Lorentzian 
and 1/f noise of several different noise levels. Figure [6] 
compares the Fano factor with leakage and without leak- 
age; in every case, leakage increases the Fano factor no- 
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FIG. 5. (color online) Fano factor [Eq. ([9|] for the non-leaky 
(perfect) integrate and fire model as a function of time for 1/f 
noise, using the same parameters as in Fig. |4] We show sim- 
ulations of 1/f noise for two different time series: T = 10^ s 
and T — 10'^ s, where T is the length of the simulation time 
window. We also show the corresponding analytical approxi- 
mations Eq. i\17\ . Note how the Fano factor due to 1// noise 
increases logarithmically until it reaches a slight saturation 
at t = T, the "maximal experimental time". This behavior 
agrees qualitatively to what is observed in measurements on 
single visual neurons in cats and insects [31j . 



ticeably (i.e. there is increased variability in the spike 
train). This happens because in the presence of leak- 
age, the neuron tends to forget past inputs that were not 
strong enough to break the threshold barrier and returns 
to its rest state even though it received a considerable 
amount of subthreshold input. This situation is dramat- 
ically different from the perfect (non-leaky) model: In 
the absence of leakage, every subthreshold stimulation 
increases the charge of the neuron, priming it for firing. 
The quantitative difference in Fano factor highlights the 
importance of leakage in neuron dynamics. For example, 
note the dramatic quantitative difference between Fano 
factor for leaky and non-leaky cases at the level of 1% of 
1/f noise {h/Io ^ 0.01). 

Our explicit numerical results of Fano factor subject to 
1/f noise shows that F(t) does not saturate at long times. 
This property is consistent with the the direct measure- 
ments of neuron Fano factor presented in Teich et. al. 
(See Fig. 6 in [31 ). Their experimental Fano factor in- 
creases well beyond one, and does not appear to reach 
a plateau at long times. Note that previous calculations 
based on Lorentzian noise [HJ [35] have showed that the 
Fano factor saturates at long times. This result is in- 
dependent of whether one uses the perfect or the leaky 
integrate and fire model. 

One can argue that with respect to a Lorentzian noise 
model, the experimental data only explicitly rules out 
cases with Tc < T, i.e. Lorentzian models with correla- 
tion time longer than the experimental window will not 
show saturation. This is true, but the Lorentzian and 
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(a)Lorentzian Noise, 7 = 1 Hz (b)l// Noise, T = 100 s 

FIG. 6. (color online) Fano factor for integrate and fire neuron with leakage and without leakage, for /i/Jo ~ 1,0.1,0.01. 
Leakage increases the Fano factor considerably, and can produce large quantitative difference at low levels of noise. 



1// model predictions for F{t) differ in other respects as 
well, which do not depend on Tc or T. 

The Fano factor F{t) for the experimental data has a 
tangent slope of 0.4 — 0.5 around t ~ 1 s in log-log plot 
(Fig. 6 of |3I])- The predicted slope of the Lorentzian 
Fano factor, however, is 1 for t ^ Tc (before saturation 
takes place) . This difference rules out the possibility that 
the experimental results can be understood in terms of 
an unsaturated Lorentzian model. 

The Fano factor under 1// noise has a tangent slope of 
approximately 0.7 around t ^ 1 s (Figs, [s] and |6]d). Note 
that the Fano factor increases logarithmically, so the ex- 
ponent depends on particular time t chosen to measure 
the slope. We repeated our simulation using 1//°'^ noise, 
and obtained results very similar to Figs. [5] and [6)3, ex- 
cept that the slope was reduced to « 0.5. This shows 
that a LIF model subject to 1//" noise can provide an 
excellent fit to experimental data of long time neuron 
dynamics, and further suggests that neuron input noise 
is better approximated by a l//-like spectrum, than a 
Lorentzian with a single characteristic correlation time. 

IV. HOW NEURONS RESPOND TO A SUDDEN 
STEP EXCITATION 



120 
^ 100 
B 80 
2 60 



Lor. y=0 Hz ■ 

1/f !!= 

Lor. y=1000 Hz r" 
White 
No Noise 1 


1 
1 








lllllllllll 


111 
llliiiiBillllllllll 



] I 'fiin m I ■ ^ 

1.48 1.5 1.52 1.54 1.56 1.58 
time (s) 



FIG. 7. (color online) Mean fire rate of single neuron in 
response to step excitation under various types of noise: 
Lorentzian noise with 7 = Hz ("static" case), Lorentzian 
noise with 7 = 1000 Hz, white noise, and 1/f noise. We also 
included the noiseless case for comparison. The presence of 
noise primes the neurons to react significantly faster. The re- 
action time is optimal (almost instantaneous) for the "static" 
case. Surprisingly, under 1/f noise the neuron response time 
is nearly optimal. 



We now analyze the impact of low frequency noise on 
the reaction time of a LIF neuron. We consider the re- 
action of a neuron subject to the following stimulus, 

I{t) = {hr,{t) + he{t ~ tstep)) ■ (18) 

For times between t — and t — igtcp: the current input 
is pure noise with subthreshold amplitude Ii;att~ tstcp, 
a superthreshold bias current /q is suddenly turned on in 
addition to the noise. In the simulations below we used 
^stop = 1-5 s, Jo = 4.3 X 10-1° A, and h = O.S/q. Fig- 
urelT] shows the mean fire rate averaged on 1 ms bins after 



100,000 time series are taken into account. The ensemble 
of time series can be thought of either an actual ensemble 
of different neurons or to a single neuron subject to sta- 
tistically similar excitations at different times. In both 
cases, we can make claims of optimality on "average" . 

In the absence of noise, the LIF neuron takes 43 ms to 
respond. The quickest response is obtained in the case 
of Lorentzian noise with 7 = Hz, that is equivalent 
to S{f) — (5(/)/(27r), the static limit discussed in sec- 
tion |II A| Note that in this case the noisy current does 
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FIG. 8. (color online) Histogram of first spike times under 
a step excitation subject to various types of noise. In the 
case of white noise, it takes 30 ms for 1% of the neurons to 
spike for the first time. The situation is markedly different 
for the static case (Lorentzian with 7 = Hz) and for 1/f 
noise. Here some neurons respond almost immediately, while 
others take a long time to spike; note the long time tail in the 
distributions for Lor. Hz and 1// noise. 



not change in time, and is equivalent to a bias current 
with amplitude picked from a Gaussian distribution. The 
response under 1/f noise lags behind the static case by 
only w 5 ms, i.e., it is nearly optimal. Both static and 
1/f noise reach their steady state much faster when com- 
pared to other types of noise. Lorentzian noise with roU- 
ofF frequency 7 — 1000 Hz is intermediate between white 
and 1/f noise. Clearly, the response time improves as 
the noise gets dominated by low frequency components. 

In order to elucidate the mechanism by which noise 
sensitizes neuron response time, we present two addi- 
tional figures. Figure [8] shows a histogram of first spike 
times. Under white noise conditions, it takes approxi- 
mately 30 ms for 1% of the neurons to spike for the first 
time. On the other hand, it takes only 8 ms for 1% of the 
1/f noise neurons to spike and a mere 2.5 ms for 1% of 
static noise (i.e. Lorentzian 7 = 0) neurons to spike for 
the first time. Neurons subject to 1// noise and the 7 = 
static noise have a much wider distribution of first spike 
times. Some neurons spike almost immediately, while 
others take a long time to spike - note the extended tail 
in the distributions. This broad distribution implies a 
greater variability in the neuronal interspike interval and 
a degradation of any information coded therein. Extreme 
low frequency noise {1/f and 7 = 0) result in a trade off 
between reliability and rapid response. 

Figure [9] sheds light on the origin of this effect, by 
plotting the distribution of neuron voltages just before 
the step stimulus is applied. Here we see why the static 
case is optimal: The distribution of neuron voltages is 
nearly flat, and extends close to threshold. Hence, when 
the stimulus is applied, a significant amount of "primed 
neurons" will reach threshold almost instantly. While the 
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FIG. 9. (color online) Neuron voltage right before the step ex- 
citation arrives. Note that for the static case (Lor. Hz) and 
for the 1/f noise case, a significant fraction of neurons have 
quite high voltage. Hence these neurons will reach threshold 
much faster, explaining how noise can "prime" neurons for a 
fast response time. 



voltage distribution for 1/f noise is not flat, it is broad 
and extends all the way to high voltages. Similar to the 
static case, the presence of a tail extending near threshold 
implies that a significant number of neurons are "primed" 
by 1/f noise; these neurons will react nearly instantly to 
the stimulus. 

We note in passing that it is possible to engineer a 
Lorentzian noise spectrum to yield a response time sim- 
ilar to that for 1/f noise via an appropriate choice of 
7 = where we recall that T is the simulation time 
window. For T = 2 s, we find that 7 = 0.5 Hz does 
the trick, in agreement with expectations based on our 
previous analysis of the distribution of interspike time in- 
tervals for different conditions (c.f. Fig. [2]): The shortest 
interspike time interval for a Lorentzian with 7 = 0.5 Hz 
is similar to that for 1/f noise. The two distributions are 
not identical. The ISI for the 1/f has a tail that extends 
to much longer times. 



V. CONCLUSIONS 

In conclusion, we presented analytical and numerical 
calculations of the perfect and the leaky integrate and 
fire neuron aimed at elucidating the impact of 1/f noise 
on single neuron dynamics. Though more difficult to 
analyze than the perfect integrate and fire model that 
is commonly used in noise and network studies, our LIF 
model is a more realistic model of a neuron and as we 
have shown, the inclusion of the "leakiness" gives rise to 
much higher disorder. 

With regard to the response of the LIF to 1// noise, 
we find a surprising dichotomy: While it degrades the 
ability to transmit information using interspike times, it 
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manages to enhance the overall response time (of an en- 
semble of neurons) to a sudden stimulus by ensuring that 
a subset of neurons are primed with a near threshold volt- 
age. 

Our explicit numerical simulations of neuronal re- 
sponse to a sudden step excitation elucidates the mecha- 
nism by which noise can enhance neuron response time. 
Neuron response times were shown to be optimal under 
static noise (noise sharply peaked at zero frequency) be- 
cause a large number of neurons are close to threshold 
just before the step stimulus arrives. This select group 
of neurons act as the alarm bells: They respond almost 
instantly to the step stimulus. Remarkably, as we demon- 
strated, the case of 1// noise is not much different; the 
distribution for neuron voltages just before stimulus also 
has a long tail extending towards threshold. These results 
allude to a possible explanation for why the brain is pop- 
ulated by an astronomically large number of neurons. It 
may well be an evolutionary adaptation designed to take 
advantage of the ambient noise to enhance the probability 
of survival. Neuron redundancy enables faster response 
times in the presence of low frequency noise, which in 
turn allows an animal to react quickly to a sudden dan- 
ger. 

However, the apparently beneficial feature noted above 
does not come without cost. The 1// noise trades off 
speed for reliability by introducing much more variability 
in the properties of the resulting spike train. We quan- 
tify this uncertainty using the Fano factor. Our analysis 
of the Fano factor reveals that in the presence of 1// 
noise, this measure of disorder increases logarithmically 
as a function of time. On a positive note, we find an 
excellent qualitative agreement between Fano factor for 
the 1// noise and the Fano factor derived from labora- 
tory results of experiments with single neurons (311 132j . 
Specifically, the latter also rises monotonically well be- 
yond one and shows no evidence for saturation. This 
agreement suggests that the neuron input noise is better 
approximated by a scale-free l//-like spectrum than the 
more commonly invoked low frequency Lorentzian spec- 
trum. 

The rate at which the Fano factor for a Lorentzian 
spectrum grows with time prior to saturation is F(t) oc t, 
independent of 7. Moreover, the Lorentzian Fano factor 
always tends to a plateau at long times [211 ES] ■ Both 
tendencies are at odds with the behavior of the experi- 
mentally measured Fano factor regardless of whether the 
Lorentzian spectrum is fine-tuned to yield a response 
time similar to that of 1// noise. We note that these 
claims can be definitively tested by repeating the exper- 
iments of Teich et al. [31j with increasingly longer exper- 
imental time windows T. 

The logarithmically rising Fano factor reflects the fact 
that the long time spike dynamics is dominated either by 
periods of extended inactivity or by periods of aggressive 
bursting. This behavior is due to the lack of ergodicity 
in 1// noise, i.e., the fact that it lacks a characteristic 
correlation time. Not surprisingly, therefore, the neuron 



dynamics in the presence of 1// noise is very different 
from that due to Lorentzian noise. As an aside, we note 
that the degree of uncertainty is also substantially greater 
in the leaky model than in the non-leaky (perfect) case. 

These conclusions are consistent with the observation 
that some neurons seem to spike in a very irregular fash- 
ion. The temporal gaps of a spike train has much larger 
information capacity and for this reason, there is consid- 
erable body of work arguing that neurons use the timing 
intervals to encode information. Loss of reliability due 
to low frequency noise, however, limits the information 
capacity of the spike trains [22] . On the other hand, re- 
search has shown that in certain cases neurons can spike 
with high degree of reproducibility [3S]. Whether the 
origin of highly reproducible spike patterns is due to ex- 
tremely low noise at the single neuron level, or due to a 
network effect that compensates for the noise remains to 
be seen. A more interesting possibility is that the various 
functional regions of the brain may have evolved differ- 
ent strategies for managing ambient noise, depending on 
function and associated information capacity demands. 
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Appendix A: Numerical simulation of Gaussian 
noise with arbitrary spectral density 

To simulate the noise used in Eq. ([2| we used a varia- 
tion of the efficient algorithm proposed by Timmers and 
Koenig [361 . Consider the time window from t — to t — 
T. Define a discrete set of N time instants tm — nA</2, 
where n = 0, 1, . . . , 2iV - 1 and At = T/N. Choosing N 
as a power of 2 allows the use of the fast Fourier trans- 
form algorithm, with significant speed up. The associ- 
ated set of "lower half" frequencies are fm = m/T, with 
m = 0,1,..., iV, and the "upper half frequencies" are 
= (2iV - m)/r for m = AT + 1, + 2, . . . , 2iV - 1. 
We are now ready to state the algorithm that generates 
individual real-valued time series rj{tm) (??(/) are their 
Fourier transforms): 

1. Set 77(/o) = 0; 

2. For each m = 1, . . . , A^ - 1, set ?? = e"''- \/s{fm), 
where is a random number in the interval [0, 1); 

3. Set 7y(/w) = y^^(/jv); 

4. Set 77(/Ar+m) equal to the complex conjugate of 
vifN-m) for all m = 1, . . . , N; 

5. Finally, take the inverse Fourier transform of fj^fm) ■ 
The resulting ri{tm) realizes an individual time se- 
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ries of the Gaussian process with noise spectrum 

Sifm). 

Figure [lO] depicts three example time series: white 
noise, 1/f noise, and Lorentzian with half- width 7 = 
10 Hz. We simulated 100,000 of these time series and 
studied their amplitude distribution and noise spec- 
tra [Eq. ^ 



Figure 
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demonstrates that the noise 
amplitudes are distributed according to a Gaussian, 
and Fig. [12] computes the ensemble average of their 
correlation function, S{t) {r]{t)r]{0)) . The latter 
have the expected forms: For white noise, S{t) = 
sin (27r7max0/(27r7,naxi), for 1/f noise, S{t) « 1 - [Cb + 
ln(7maxt)]/ln(7max/7min) {C E = 0.5772 is the Euler- 
Mascheroni constant), and S{t) — e^^* for Lorentzian 
noise. 



Appendix B: Exact calculation of the interspike 
interval histogram for the LIF model subject to 
static noise 

The case of static noise S{f) = (5(/)/27r (a Lorentzian 
with 7 — 0) is particularly simple because the stochas- 
tic process rj{t) randomizing the current Eq. ^ does not 
change in time. Each 77 is picked from a Gaussian distri- 
bution at t = 0, 



In the presence of leakage, the voltage is obtained by 
solving Eq. ([T]) under a constant current I = Iq + Iirj, 



V{t) = R{Io + hv) (1 



e «c 



(B2) 



The interspike time interval I will be given by the time 
it takes for this voltage to reach Vth, leading to 



I = Tr-RCln 1 - 



Vtb/R 



(B3) 



where is the refractory time period. The ISI distribu- 
tion can now be computed from the expression 



P{1) 



1 



dr/p{ri)d 



I - [Tr.- RCln 1 



'»? J-00 

where the normalization factor 



[Vtb/{RIi)-Io/h] 



lo + hv 
(B4) 



(B5) 



ensuring that rj is strong enough to "click" the delta func- 
tion. After some algebra we obtain the following exact 
result: 



(Bl) 



As a result, the quasi-static method for calculating the 
ISI distribution developed in Refs. [H] and jM] becomes 
exact. As we show below, this allows us to compute the 
ISI distribution exactly even in the presence of leakage 
{R < 00). 



P{1) = 



1 lo' 

2 Wi . 



RC 1-0 RC 



2-^1^'^ [i?c(l-e-^) 



This expression is plotted in Fig. |2] 



[1] C. Koch, Biophysics of Computation: Information pro- [11 
cessing tn single neurons (Oxford University Press, Ox- 
ford, UK, 1999). [12 

[2] M.B. Weissman, Rev. Mod. Phys.60, 537 (1988). 

[3] Sh. Kogan, Electronic Noise and Fluctuations in Solids [13 
(Cambridge University Press, Cambridge, 1996). 

[4] R.L. Badzey and P. Mohanty Nature 437, 995 (2005); [14 
A.R. Bulsara, ihid. 437, 962 (2005). [15 

[5] R. Almog, S. Zaitsev, O. Shtempluck, and E. Buks, Appl. 

Phys. Lett. 90, 013508 (2007). [16 

[6] K. Murah, S. Sinha, W.L. Ditto, and A.R. Bulsara, Phys. 

Rev. Lett. 102, 104101 (2009). [17 

[7] D.N. Guerra, T. Dunn, and P. Mohanty, Nano Letters 9, 

3096 (2009). [18 

[8] D.N. Guerra, A.R. Bulsara, W.L. Ditto, S. Sinha, K. 

Murah and P. Mohanty Nano Letters 10, 1168 (2010). [19 

[9] J. Zamora-Munt and C. MasoUer, Optics Express 18, 

16418 (2010). [20 

10] L. Worschech, F. Hartmann, T.Y. Kim, S. Hofling, M. 

Kamp, A. Forchel, J. Ahopelto, I. Neri, A. Dari, and L. [21 
Gammaitoni, Appl. Phys. Lett. 96, 042112 (2010). 



P.I. Fierens, S.A. Ibaiies, R.P.J. Perazzoa, G.A. Patter- 
sonb and D.F. Grosza, Phys. Lett. A 374, 2207 (2010). 
G.L. Gerstein and B. Mandelbrot, Biophys. J. 4, 41 
(1964). 

A. L. Hodgkin and A. F. Huxley J. Physiol. 117, 500 
(1952). 

L. Lapicque, J. Physiol. Paris 9, 620 (1907). 

A. Longtin, A. Bulsara, and F. Moss, Phys. Rev. Lett. 67, 

656 (1991). 

J.K. Douglass, L. Wilkens, E. Pantazelou, and F. Moss, 
Nature 365, 337 (1993). 

Z.F. Mainen and T.J. Sejnowski, Science 268, 1503 
(1995). 

D. Nozaki, D.J. Mar, P. Grigg and J.J. CoUins, Phys. 
Rev. Lett. 82, 2402 (1999). 

D. Nozaki, J.J. Collins, and Y. Yamamoto, Phys. Rev. 
E60, 4637 (1999). 

R. Soma, D. Nozaki, S. Kwak, and Y. Yamamoto, Phys. 
Rev. Lett. 91, 078101 (2003). 

J.W. Middleton, M.J. Chacron, B. Lindner, and A. 
Longtin, Phys. Rev. E68, 021920 (2003). 



12 



T3 
-I—' 

■q. 
E 
< 




0.05 
time (s) 

(a)White Noise 



0.1 



T3 
-I—' 

■q. 
E 
< 




0.05 
time (s) 

(b)l// Noise 



0.1 



■D 



0.05 
time (s) 

(c)Lorentzian Noise, 7 = 10 Hz 



0.1 



FIG. 10. (color online) Example time Series for three different kinds of noise spectra. 
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FIG. 11. (color online) Histogram for the amplitudes of r]{t) for 100,000 time series, demonstrating that ri{t) is Gaussian 
distributed. 
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FIG. 12. (color online) Ensemble average correlation function, S{t) — (?7(i)»7(0)), calculated as an arithmetic average of several 
time series generated by the method. 



